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Abstract 

Simple nonlinear climate models yield a precession index-like term in the 
temperature. Despite its importance in the geologic record, the precession index e sinrn, 

where e is the Earth’s orbital eccentricity and mis the Sun’s perigee in the geocentric frame, 

is not present in the insolation at the top of the atmosphere. Hence there is no one-for-one 
mapping of 23,000 and 19,000 year periodicities from the insolation to the paleoclimate 
record; a nonlinear climate model is needed to produce these periods. Two such models, a 
grey body and an energy balance climate model with an added quadratic term, produce <? sin 

m terms in temperature. These terms, which without feedback mechanisms achieve 

extreme values of about ±0.48 K for the grey body and ±0.64 K for the energy balance 

model, simultaneously cool one hemisphere while they warm the other. Moreover, they 
produce long-term cooling in the northern hemisphere when the Sun’s perigee is near 
northern solstice and long-term warming in the northern hemisphere when the perigee is 
near southern solstice. Thus this seemingly paradoxical mechanism works against the 
standard model which requires cool northern summers (Sun far from Earth in northern 
summer) to build up northern ice sheets, so that if the standard model is correct it may be 
more efficient than previously thought. Alternatively, the new mechanism could possibly 
be dominant and indicate southern hemisphere control of the northern ice sheets, wherein 
the southern oceans undergo a long-term cooling when the Sun is close to the Earth during 
southern summer. The cold water eventually flows north, cooling the northern 

hemisphere. This might explain why the northern oceans lag the southern ones when it 
comes to orbital forcing. 
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1. Introduction 


An important component of Milankovitch’s [1941] astronomical theory of climate 
change is the precession index. The precession index, along with the Earth’s obliquity and 
orbital eccentricity, are believed to be the major controling factors of climate change in the 
last few million years [e.g., Hays et al, 1976; Berger et ai, 1984; Hinnov, 1999]. The 

precession index is e sin <y, where e is the Earth’s orbital eccentricity and (assuming a 

geocentric point of view) 0) is the argument of perigee of the Sun’s orbit about the Earth. 

The precession index spectrum has a major peak at 23 kyr and a smaller one at 19 kyr (1 
kyr = 1000 years). 

It is an interesting fact that the equation for the insolation at the top of the 
atmosphere contains no terms which look like e sin co ( Rubincam [1994; 1996]; Hinnov , 

1999; and Bruce G. Bills, private communication, 1994). Therefore the equation contains 
no long-period 23 kyr and 19 kyr terms. Because of its undoubted importance in the 
paleoclimate record, the existence of the precession index must be due to the Earth’s 
nonlinear response to the insolation. 

The standard explanation for the importance of the precession index, namely that 
cool northern summers are required for the growth of ice sheets [e.g., Milankovitch, 1941, 
pp. 435-435] is in fact such a nonlinear model. Assuming the Earth responds to mainly to 
northern summertime cooling is to manipulate the insolation into giving 23 kyr and 19 kyr 
frequencies in ice volume [Rubincam, 1996]. However, the standard “model” given by 
Milankovitch [1941] is not really a quantitative model; instead he merely correlates the 
snow line with summer insolation and calls the result “proof.” This procedure is 

questionable, given all of the factors which might contribute to the position of the snowline 
besides just insolation. 
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Below I demonstrate below two quantitative nonlinear models which produce e sin 
ft) in temperature from the e sin a>free insolation. The first is the simplest of all “climate” 
models: a grey body. In this model insolation is proportional to T 4 , where Tis surface 
temperature. I write temperature as T = T 0 + AT and expand T 4 to order (AT) 2 . The other 
is an energy balance climate model with added T 2 , T\ and T 4 terms; these are also 
expanded to order (AT ) 2 . In both these nonlinear models cross-product of terms in the 

insolation produce AT pi « -e sin 0)P,(sin <f>) in the surface temperature, where P/sin (j>) is 
the Legendre polynomial of degree 1, 0 is latitude, and the subscript “pi” on AT stands for 
precession index”. The magnitude of the temperature change can reach about ±0.48 K for 
the grey body and about ±0.64 K for the energy balance model when the eccentricity e 
reaches its maximum value of 0.06 and co= 90°, 270°. 

Because P/sin <p) = sin <p, which changes sign when the equator is crossed, these 

models indicate simultaneous wanning in one hemisphere and cooling in the other. 
Moreover, both models yield a precession index with a seemingly paradoxical sign: long- 
term cooling in the northern hemisphere when the Sun is near perihelion during northern 
summer and a long-term warming in the north when the Sun is near aphelion during 

northern summer. This effect thus works in opposition to the usual explanation of e sin co , 

wherein cool summers (Sun far from Earth during northern summer) are required for the 
snow linger in order to build up into ice sheets. In this case the new mechanism warms the 
ground during the cool northern summers, making it harder for the standard model to 
work. 

Alternatively, the effect presented here may dominate the standard mechanism. The 
new effect might argue for southern hemisphere control of the northern ice sheets, the idea 
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being that the southern oceans cool when the Sun is near the Earth during southern 
summer; the cold water flows north, eventually cooling the northern hemisphere enough to 
build up ice sheets. This scenario might explain why the northern oceans lag the southern 
in the orbital forcing. 


2. Insolation 


The insolation is given by 


F s = QS 

s 


where 



with F° being 1371 W m' 2 when the reference distance r 0 is r 0 = a = 1 AU [Hickey et al., 
1988], where a is the semimajor axis of the Earth’s orbit, and 


5 - 


1=0 m = 0 
l +« 


9 S q (e) 

p=0 q=-o o 

\cn<iY~ meven 

*[sin 1 _ modd [(^—2/7)0 + (£-2p+q)M + 


( 1 ) 
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[Rubincam, 1994]. Here e is orbital eccentricity, mis the argument of perigee, M is the 
mean anomaly of the Sun’s orbit about the Earth (I assume a geocentric point of view), and 
Q is the position of the line of nodes of the Sun’s orbit measured from the equinox; Q - 0. 

The angle e is the Earth’s obliquity (currently about 23.44°), and <p is latitude and X is east 
longitude on the Earth, while t\ is the hour angle. The PJsin <f>) are the associated 

Legendre polynomials of degree /and order m . The FJe) are the inclination functions 

from celestial mechanics; Rubincam [1994, Table 2] gives them for degree 1, while Kaula 
[1966] and Caputo [1967] list them for degrees 2 through 4. The W r2pq (e) are not the 
eccentricity functions from celestial mechanics; rather, they are special eccentricity 
functions associated solely with the insolation and are tabulated in Rubincam [1994]. 

(There is a typographical error in Table 3 of that paper; the entry for /- 2p = ±1, q = ±1 

should read 2e -3e 3 /2 instead of 2e - 3e 4 /2. Also, the (2'- 2s)\ factor in the equation 
following (7) in the same paper should read (2/- 2 j)!) 

The zonal insolation can be found from (1) by setting m = 0. The zonal insolation 
given by equation (8) in Rubincam [1994] contains some errors. The corrected equation 
for the zonal insolation at the top of the atmosphere is 


-fc*. | Ln 
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^o(sin <t>) • [l + + 2e cos M + -e 2 cos2A/l 

2 2 J 

+ 2 Pj (sin <f>) 

* - ~^) s i° £ sin ( co + M ) + 2e sin e sin (co + 2 M) 

+ ~S Sin C Sin + ~^e 2 sin e sin (co + 3A/)J 

P (sin <j>) 

# [( 1 + y)(|sin 2 £-|) - ^(1 - ^ 2 ) sin 2 e cos (2co + 2M) 
3 1 3 

+ 2e(-sin 2 £--)cosM + -^e sin 2 e cos (2co + M) 

9 r <5 * 

— — ^sin ecos(2co + 3M) + —e 2 (—sin 2 £ — )cos2 M 

^ 2 4 2 

- -^-e 2 sin 2 £ cos (2a> + 4Af)J 


+ ... 


The P/sin <j>) are the usual Legendre polynomials 
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P x (sin 0) = sin 0 , 


(sin 0 ) 



^3 (sin 0) = 



—sin 0 , 
2 v 


etc. Being just the zonal (longitude-independent) insolation, the diurnal terms are not given 
in (2). 


The important thing to notice about (1) and (2) is that they do not contain a term of 
the form 


S - ... + coefficient x e sin co + ... 

In other words, there is no term in the insolation which looks like the precession index. 
Hence there is no term in the insolation which has periods of 23 kyr and 19 kyr. Bruce G. 
Bills, in an independent analysis (unpublished), arrived at the same conclusion. This lack 
of a precession index was in fact demonstrated early on by Humphreys [1964, pp. 85-87] 
in his classic book on meteorology. 

The only thing which even remotely looks like the precession index is the 2e sin e 

sin (a>+ 2M) term which multiplies P, (sin 0). However, this term has a high frequency 

due to the presence of the 2 M. Its period is about half a year, not 19 kyr or 23 kyr. 

This result has caused some confusion [Berger, 1996; Rubincam , 1996] because 23 
kyr and 19 kyr periods are strongly indicated in the paleoclimate records [e.g.. Hays et al . , 
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1976, Berger et al . , 1984]. If the Milankovitch theory of climate is correct, how can there 
be an e sin (o signal in the paleoclimate records when it does not exist in the insolation? 

The answer is that if e sin co is important for climate, it must be due to the way the 

Earth responds to the insolation [Rubincam, 1994, 1996]. In other words, the Earth’s 
climate system does something nonlinear to the astronomical signal, thereby manufacturing 
the 23 kyr and 19 kyr periods. Rubincam [1994, p. 201] produced a model which in fact 
gave a precession index-like term in the radiation reaching the ground. That radiation is the 
insolation at the the top of the atmosphere multiplied by the Earth’s coalbedo. Short-period 

terms in the albedo can multiply short-period terms in the insolation, eliciting e sin co. But 
this precession index” was extremely weak. Below I exhibit two simple models which 
which produce e sin to in the Earth’s temperature and have extreme values of about ± 0.48 

K and ± 0.64 K when the Earth’s orbital eccentricity is at its maximum value of e = 0.06 
and co = 90°, 270°. 

3. Grey body ' 

The first is the simplest of all “climate” models: the Earth as a grey body with 
albedo A and emissivity e a . In this case 


£„ a T 4 = (l-A)QS . 


( 3 ) 


where <7 is the Stefan-Boltzmann constant. Writing 
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T = T 0 + AT 

where T is the absolute temperature yields 


10 

(4) 


o 4 


1 + 4 


AT) 

+ 61 




+ 4| 


t) 

\*0 ) 


AT 


= (1~A)QS = (1 -A)Q + (1 - 


(5) 


where T 0 is a constant and 


AS = 2P 0 (sin 0) • £ e cos M~\ 

+ 2 P, (sin 0) • [(sin e sin (a) + M) + 2e sin £ sin (o) + 2M)J 
+ ^ 2 (sin^)fjsin 2 e-i-')*[l + 2eco sAf] 


( 6 ) 


to the first power in e in (2). Setting 


= (1 -A)Q 
x = (A T/T 0 ) 


y = (1 - A) QAS 
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in (5) makes that equation become 


y = 4x + 6x 2 + 4x 3 + x*. 


This equation has the form 


y = a x x + Oj x 2 + Oj x 3 + a 4 x A + ... 


Assuming 


AT « T 0 


(8) can be solved for * using reversion of series [Selby, 1974, p. 470]: 


x 


Ay + Ay 2 


+ A y 3 + Ay* + — 


where 


A, = 1/fl,, Aj = -flj/o, 3 , 


Hence the solution of (5) becomes 


AT 


(l-A)QM 


_3_ 

32 


(1 ~ AQAS 

£ ss 


n2 
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Squaring AS as given by (6) in the above equation yields cross-product terms of the form 


e sm e sin co 


^ p ,(s in*) - P i (sin <f>)P 2 (sin <f>) 


to first power in sin e and are mainly from the e cos M multiplying sin e sin (co + M) P/sin 
<(>). This is easily seen from the identity 


sin a cos/J = sinja + p) 

2 2 

where a= co + M and /3 = M, giving a low-frequency sin to and a high-frequency sin (co + 
2M). The high-frequency terms will be ignored here. Using 

2 3 

/ , 1 (sin0)/> 2 (sin0) = -/>(sin <j>) + -/> 3 (sin0) (10) 

ultimately gives 


AT 

T 0 


9_ (1-A)Q 

32 L 


• esme sin co 


-Pi (sin <p) - j P 3 (sin 0)J 


+ other terms 
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Hence the “precession index” AT pi in temperature is 


AT- = - 


P l 


9T 0 \ 

'(1 ~A)Q 

32 



12 


• e sin e sin co 


P i (sin 0) 


— P j(sin 0), 


( 11 ) 


for a grey body for zonal insolation only. For the Earth the variations in e are small and 
can be considered constant to a first approximation, so that this term behaves like e sin co. 
Using <t= 5.67 x 10' 8 W m' 2 K -4 , together with some typical values A = 0.3, r 0 = a = 1 

AU, F s ° = 1371 W m‘ 2 , e = 23.44°, plus setting e u = 0.89 with the emissivity e a being 

chosen to reproduce the observed average surface temperature of the Earth T 0 = 288. 1 K 
[North et al., 1981, Table 1] give 


^pi ~ -5.6»esin<2)|^ P,(sin0) - -jPjCsin^J K (12) 

Ate = 0.06, the maximum obliquity of the Earth’s orbit, and at latitude 0 = 43.1°, this 
reaches a minimum value of AT pi = - 0.27 K when co= 90° (perihelion at northern 
solstice). At the same time at 0 = - 43. 1 ° it reaches a maximum of AT pi = + 0.27 K. The 

curve marked “GREY (Z)” in Figure 1 shows the temperature as a function of latitude for e 
= 0.06 and co= 90° for a grey body using only the zonal insolation. 
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4. Grey body by iteration 

It will instructive to solve for AT pi again through an iterative process, since it will 

be the method used to solve the energy balance equation below. Reversion of series cannot 
be used to solve a differential equation like the energy balance equation. 

In this case the first-iteration temperature will be written 

T = T 0 + AT . 


To first order (5) is simply 




( at) 

— = (l-A)QAS 

\*0 J 


after using (7) to alias out the constant terms. To second-order 


T = T, + AT + ST 


so that (5) will now be appro xim ately 


^ss G T 0 


(AT +3^ J AT^ 2 


+ 6 




V r o J 


= (1 -A)QAS 


(13) 


The (AT) 2 part will yield 
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(A7) 2 


3 _ 

16 


' (l-A)g T 


• e sin £ sin co 




P i (sin 0) 


2 


(sin <t>) 


plus other terms. Since the terms in the above equation not exist on the right side of (12), 
they must be cancelled by similar terms in ST ; i.e.. 


4 



+ 



= 0 


for the e sin fit) terms, thus giving 


ST . = 

p< 


9T 0 


32 


(1- A)g T 


• e sine sin fit) (sin - ^-P 3 (sin 


which is the same as (1 1). This iterative approach will be used in the next section. 

The above equations for the grey body uses only the zonal insolation (2). 
However, because the grey body has no thermal inertia, the diurnal terms in (1) will be of 
the same order of magnitude as the zonal terms. Including the /= 1,2 m = 1 terms in the 
above approach yields the curve marked “GREY (Z + D)” in Figure 1. The diurnal terms 
increase ther size of the effect to extremes of ±0.48 K. 
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5. Energy balance climate model 


The zonal nonlinear energy balance climate model is [e.g.. North et al, 1981, 

1983]: 



2 .dr 


. ty 


dll 


+ A) + X fi . r = (1-A)es 

1=1 


(14) 


where \1 = j/n <j> and the B t T‘ terms for i > 2 are the added nonlinear terms, and where now 

Tis temperature in degrees Centigrade instead of the absolute temperature. Qnly the zonal 
terms will be needed, since the short-period diurnal terms will be damped out by the heat 
storage term which was not present for the grey body. 

The first term CdT/dt in the equation above represents heat storage, while the 

second term represents zonal diffusion as produced by winds and ocean currents and has a 

diffusion constant D. The infrared radiation leaving the Earth in the above equation is 
given by 


/ = A) + ByT + B 2 T 2 + B 3 T 3 + B,T 4 (1: 

The values for A 0 , and the B : can be estimated from the data of Graves et al. [1993] (see 
Figure 2); these values are A 0 = 195.0, B, =1.4158, B 2 = 0.02289, B 3 = 0.001 148, and 
B 4 = 0.00002089. These numbers do not come from a least-squares fit; rather, they were 
chosen to bisect the envelope of the data. 

The linearized version of (14) will be solved for first. Once again writing 
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T = T 0 + AT 

gives by (15) 

I = Ho + 

1=1 

where 

H o = + B X T 0 + b 2 t 2 0 + b 3 t 3 + b 4 t 4 0 

H\ = B^ + 2 B 2 T 0 + 3 B 3 Tl + 4 B 4 T 3 0 
H 2 = B 2 + 3B 3 T 0 + 6 B 4 Tl 

H 3 = B 3 + 4B 4 T 0 

h 4 = B 4 

and now T 0 = 14.9° C. [North et al., 1981, Table 1], with H 0 = 226.005, H, = 3.1389, H 2 

= 0.1020, H 3 = 0.00235, and H 4 = 0.00002089. The linearized, first-order equation is 
then 
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r dO^ + AT) 

dt 


D(l-n 2 ) 


u d(T 0 + AT) 


dp 


dp 


+ H 0 + H X AT = (1-A)0[1 + AS] 


Setting 


H 0 = (1 - A)Q 


allows the albedo A to be solved for, giving A = 0.344, in reasonable agreement with the 
observed value of 0.30 [Stephens et ah, 1981], and leaving 


r jK AD 

dt 


d\ 


D(l-p 2 ) 


dp 


dp 


+ H x AT = (1-A)QAS . 


(16) 


In this equation AT is assumed to have the same form as the zonal AS as given by (6): 

AT = r^2ecos (M - y/) 

+ 2 P x (sin 0) sin e sin (ta + M - y/) 

( 17 ) 

- — e P (sin 0) cos ( M - yr)J 
+ | h ^(sin 0)[|sin 2 e - i] 
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where the observed values are r, = 19.5 °C., r 2 = 58.7 °C„ and ^=31.5° [North etal , 

1981, Table 1]. (The amplitude r, and lag angle i/z'were actually found by North et al. only 

for the 2 P/sin (p) sin eterm; but they are assumed here to apply to the terms proportional to 
€ as well, since they have very nearly the same frequency.) 

The diffusion constant D can be found by equating the time-invariant parts of the 
P 2 (sin <f>) coefficients in (16), so that 


D 


(1- A)g-i7, t, 
6t 2 


0.1185 


(18) 


Finding C in (16) from the time-variable P,(sin <f>) terms is more complicated. Writing the 
Pj(sin <p) part of AT in the usual complex notation 


r [2 /’(sin <f > ) sin e e 2lt “ ] 
and substituting in (16) yields 

2^/Cf, + 2Df, + H x f, = (1 -A)Q 

where &> + M ~2n (with time 1 measured in years), so that 

t, = 


where 
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d~A)Q 

V(2D+//,) 2 +(2^C) 2 


and 


2D + 

Thus there are actually two equations involving C. Solving each equation for C using the 
observed values of H„ % and yrand the value for D found above yields inconsistent 
results: in one case C = 1.76 and C = 0.33 in the other. Fortunately the value for C is not 
needed in what follows below, and I will simply adopt the observed values of t ; and yf in 
what comes next. 

Proceeding to the second-order solution, once again I write 
T = T 0 + AT + ST 

where <5Tis the second-order part of the temperature. Substituting this in (16) and eliding 


the first-order solution leaves 
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r ® 

dt 


d\ 


D( l-y 2 ) 


dtfT) 

dll 


dll 


+ ^ 5T + H 2 (AT) 2 
= (1 - A) QAS 


(19) 


as the second-order equation. Now the H 2 (AT) 2 term produces cross-product terms which 
can be written 


3H 2 r 2 (e sine sin co) [P x (sin <p) - ~ P 3 (sin <}>)] (20 

plus other terms. Since terms of the sort (20) do not appear on the right side of (19), they 
must be cancelled by similar terms in <5T in the other parts of (19), just as in solving the 

grey body by iteration. The expression (20) above contains the long periods of e sin co, so 

that CdfSrydT changes slowly in (19) and can be neglected. To see this, assume 8 T can 

be written in the form e 2 ^, where/is the frequency associated with the 23 kyr period. 

Differentiating 8T with respect to time brings the/out in front of the term. Since /is small, 

the time derivative is small and can be ignored. Hence the precession index part of the 
temperature can be solved for from 
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. Bn _ 

dn 


+ H x 8T pi 


■ 3 H 2 Tf (e sin e sin co) [P | (sin 0) - j/^sin^)] 


producing 


ST pi = - e sin e sin co 


'mi 


2 D+H, 


r ]{/>,(« 


(sin <f > ) 




Substituting numerical values in (20) yields 


8T pi = - 1 3. 7 e sin co P x (sin 0) 


( 22 ) 


+ 5.1 esinmPjtsin^) °C. 


For e equal to the maximum value 0.06, 8T pi achieves extreme values of ST pi = ± 0.64° 

when co = 90°, 270° and (p = ± 48.5°. Thus this effect is larger than for the grey body 
(12). The results for the energy balance climate model are labeled “EBM” in Figure 1 


6. Comparison to obliquity changes 
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This result for the energy balance model can be compared to the AT obl expected 
from obliquity changes. From (17) the second harmonic in temperature is 

|T 2 F(sin *)[| sin 2 *-!] 

This will give changes in temperature as the obliquity changes with the 41 kyr period of 
magnitude 


sin e cos 


= !r 2 f|-sin ecos* |A* 

s. ^ 


where Ae is the change in obliquity. Because Ae~± l°over the obliquity cycle, 


A T m = 40.2 Ae « ±0.7 K, 

so that the precession index is about 90% the size of the obliquity variation for the energy 

balance model. These calculations of course assume no feedback mechanisms, such as ice- 
albedo feedback. 

7. Diffusion constant D 

The value derived here for the diffusion constant, D = 0.1 185, is much lower than 
the 0.649 found by North et al. [1981, p. 96] for their energy balance model; thus the 
present nonlinear model indicates that the overall heat transport is countered by other 
factors when averaged around the whole globe. The grey body gives much the same result; 
the linearized grey body model yields from (5) and (6) 
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3. 2 n 

4 Sln e _ 2^1 ^ 2 ( sin 0) = - 24.7 P 2 (sin 0) K 

for the second degree harmonic temperature variation. There will also be a second-order 
term from squaring the P,(sin <p) sin (co + M) in (6). That yields 


AT = - 


5 

(1 ~A)Q 

16 



ST 


1 

fa - a)q 

00 

1 

Q 


12 


sin 2 eR 2 (sin <p) 


-2.7 (sin 0) K 


Adding these together gives the grey body’s coefficient as - 27.4 K, which is not far from 
the observed value of - 28.0 K [North et al, 1981, p. 100.]. Since the grey body has no 
latitudinal heat transport at all, this lends support to the idea that diffusion is balanced by 
something not covered by this simple model, yielding a small D. 

On the other hand, the larger North et al. [1981] value for D indicates that latitudinal 
heat flow is substantial. And clearly diffusion has importance for the real Earth: thanks to 
the Gulf Stream London has an temperate climate, even though it is farther north than 
Winnipeg. Why are the results so different? 

The resolution to the problem presumably lies in continentality, albedo, and the way 
the infrared radiation is handled in the linear energy balance models. While the Gulf 
Stream may heat Europe, continental interiors get cold; also, the albedo increases with 
latitude [e.g., Stephens et al., 1981]. These factors, which are not included in the 
nonlinear models presented here, apparently average out the effect of diffusion, so that 
overall the Earth responds to insolation much like a grey body with a uniform albedo. The 
discepancy with the linear energy balance model may come from the linear model’s fitting 
an inherently curved set of data to the form A 0 + BT, resulting in a large diffusion 
coefficient D, thus giving diffusion a spurious importance. 
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8. Discussion 

The precession index in temperature reveals itself as the product of short-period 
terms in these T models. This shown in Figure 3 in a highly schematic diagram of the 
insolation spectrum and the temperature spectrum. The Earth’s nonlinear climate system 

creates the precession index line from the insolation. It is chiefly the e cOs M Pj[sin <p) 

term in S multiplying the strong annual term sin e sin (co + M) P/sin tp) which produces e 

sin coPj(sin (p), giving the long periods of 23 kyr and 19 kyr. It is extremely interesting to 
note that as far as long periods are concerned, the only cross-product terms of any 
significance are these e sin co terms: there are no e sin 2co or pure e terms, for instance. 

The magnitude of the effect is fairly large: for the energy balance model it is ±0.64° 
C. at maximum eccentricity and without feedback. For comparison, the difference in the 
global temperature between the present and the last glacial maximum was about -4° C. 
[Crowley, 1983, p. 868]. 

The energy balance model gives a larger temperature change than for the grey body. 
This is because the data relating infrared radiation to temperature (Figure 2) shows a greater 
curvature than for the standard T 4 behavior of a black body. 

The sign of the new mechanism appears paradoxical: when cd= 90° in (12) and 
(22), the Sun is close to the Earth at northern solstice, producing a long-term cooling in the 
northern hemisphere with simultaneous warming in the southern. Similarly, when a) = 

270°, the Sun is far from the northern hemisphere during northern summer and this effect 
produces a long-term warming there, while at the same time the southern hemisphere cools. 



Rubincam 17 December 1999 


26 


This is counter to intuition, which says that when the Sun is close to the northern 
hemisphere in northern summer that hemisphere ought to warm, not cool. However, it is 

the short period term sin e sin (co+M) P/sin (p) which does the warming in accordance 

with intuition, as can be seen from (2). The peculiar sign of the mechanism presented here 
is in fact correct and necessary to achieve energy balance in the nonlinear models. 

The sign is also counter to that of the standard explanation of why the precession 
index is important. The standard model relies on the short period terms: it calls for cool 

northern summers when co is near 270° and the Sun is far from the Earth, so that snow 

lingers through the summers, ultimately building up into an ice sheet [e.g., Milankovitch, 
1941, pp. 435-436; Pisias and Imbrie, 1986, p. 45]. By warming the northern cool 
summers and heating the cool winters, the effect found here presumably makes it harder for 
the standard model to operate. If the standard explanation is responsible for the waxing 
and waning of ice sheets at the 23 kyr and 19 kyr periods, then it must be more efficient 
than previously thought to overcome the present mechanism with its opposing sign. 

It may be the present mechanism which in fact predominates. While Imbrie et al. 
[1988] argue for northern hemisphere control of the northern ice sheets, the effect 
presented here perhaps argues for southern hemisphere control of the northern ice sheets. 

The idea is that during the times when mis near 270°, the southern hemisphere undergoes a 

long-term cooling, especially the Antarctic ice cap and the southern oceans. The cold water 

makes its way north, cooling the whole Earth and eventually producing the northern ice 
sheets. 

This could perhaps explain the phase shift seen in the sea surface temperatures for 
the precession index, as discussed by Imbrie et al [1988, pp. 144-148]. They note that the 
northern oceans significantly lag the southern oceans with respect to the orbital forcing. 

This is what would be expected from the southern oceans cooling first, and then the 
northern as the cold waters spread north. Because the tum-over time of the oceans is 
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thousands of years, the northern hemisphere should lag the southern by this amount. This 
is what is observed. There is some support for this point of view in the obliquity forcing 
as well. Here also the northern oceans lag the southern in the obliquity forcing [Imbrie et 
al., 1988, pp. 147-148], again arguing for southern ocean control of the ice ages. 

The slow ocean currents are of course not present in the model, nor are continents 
and ice-albedo feedback [ North et al., 1981, 1983; Graves et al., 1993; Short et al., 1991], 
The continents and feedback presumably amplify the effect locally. Cross-product diurnal 
terms are not present either, which could possibly give rise to significant long-period terms 
in a model better than the energy balance model. These all represent avenues for future 
research. 

The <57\ found here may not be “the” precession index. Some other nonlinear 
mechanism rather than the one discussed above may the reason for the importance of e sin 
ox The point to be made here is that there is no one-for-one mapping of the 23 kyr and 19 

kyr cycles from insolation to the paleoclimate record. The 23 kyr and 19 kyr cycles do not 
exist in the insolation. They must be manufactured from the astronomical signal by the 
Earth s nonlinear climate system. The T 2 models are one way of doing it. There may be 
others. 
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Figure 1 

AT pi for the grey body with the zonal insolation alone (dotted line), the grey body 

with the zonal and diurnal terms (dashed line), and ST p ,f or the energy balance model (solid 

line) as a function of latitude q> for e = 0.06 and 0) = 90°. These values give the maximum 
changes in temperature. 

Figure 2 

The outgoing infrared radiation as a function of temperature, based on Figure 1 of 
Graves et al. [1993], The curve, which is a polynomial in T, is chosen to bisect the 

envelope. It is the curvature of the data which gives rise to e sin co in temperature. 

Figure 3 


Schematic spectrum of insolation and temperature. This diagram shows how the 
Earth’s nonlinear climate system manufactures a new spectral line in temperature from 
existing lines in insolation. The precssion line is shown only as a single peak, and many 
other spectral lines are omitted for clarity. 
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